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ABSTRACT 

This paper is a study of a highly decomposable algorithm use- 
ful for the parallel generation of a contour surface display. The 
core component of this algorithm is a two-dimensional contouring 
algorithm that operates on a single, 2x2 subgrid of a larger grid. 
A model for the operations used to generate the contour lines for a 
single subgrid is developed. The inadequacies of the currently pub- 
lished algorithms, with respect to contour line generation for a 
subgrid, are pointed out in a brief review of the available litera- 
ture. A data structure, the contouring tree, is introduced as the 
basis of a new algorithm for generating the contour lines for the 
subgrid. The construction of the contouring tree, and the com- 
pleteness of an algorithm based upon the contouring tree, within 
the constraints of the contouring model, are shown. 

Categories and Subject Descriptors: 1.3.3 [ Picture/Image Genera- 
tion ]: surface visualization: 1.3.5 [ Computational Geometry and 
Object Modeling ]: data structures, planar contours, surface 
approximation, surface generation, surface representation, sur- 
faces; 1.3.7 [ Three-Dimensional Graphics and Realism ]: line draw- 
ings, line generation algorithms, surface plotting, surface visualiza- 
tion, surfaces; 

General Terms; contouring, contouring tree, contour surface 
display generation; 



1. Introduction 

Contour surface display generation is one of the most frequently used 
graphics algorithms [Zyda, 1984], [Zyda.1983], [Zyda, 1982], [Zyda, 1981], 

[Barry, 1979], [Faber, 1979], and [Wright, 1979]. A contour surface is a visual 
display that represents all points in a particular region of three-space <x,y,z> 
which satisfy the relation f(<x,y,z>)=k, where k is a constant known as the con- 
tour level. The function f represents a physical quantity which is defined over 
the three-dimensional volume of interest. The visual display created by this 
algorithm is the collection of lines that belong to the intersection of both the set 
of points that satisfy the relation f(<x,y,z>)=k, and a set of regularly spaced 
parallel planes that pass through the region of three-space for which the 
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Figure 1 

Contour Surface Display Generated from a Hydrogen Atom 
Wavefunction Squared (3dxy orbital) 
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relation is defined. 

For this study, the function f is approximated by a discrete, three- 
dimensional grid created by sampling that function over the volume of interest. 
The three-dimensional grid contains a value at each of its defined points that 
corresponds to the physical quantity obtained from the function, i.e. the value 
associated with point (xg.yQ.ZQ) is Vq, where f(xQ,yQ,ZQ)=VQ. In order to minimize 
confusion, we will specify the value at a particular grid point (x,y,z) by a(x,y,z), 
and will specify the value at a particular point (x,y,z) of the function by f(x,y,z). 

The visual display of the contour surface is created from this three- 
dimensional grid by taking two-dimensional slices of the grid, and constructing 
the two-dimensional, planar contours for each slice at the designated contour 
level. A slice of a three-dimensional grid is a planar, orthogonal, t\vo- 
dimensional grid assigned a constant coordinate in three-space, i.e. an x-y slice 
of a(<x,y,z>) corresponds notationally to a(<x,y>) for a particular z coordinate. 
The two-dimensional, planar contours created are the lines that satisfy the rela- 
tion a(<x,y,z>)=k for a particular planar coordinate, either x, y, or z, where 
again k is the constant contour level. If we contour all x-y slices of the three- 
dimensional grid at contour level k, we will have a stack of parallel contours 
approximating the contour surface, each planar set of contours corresponding 
to a particular z coordinate. If we contour all x-z slices of the three dimensional 
grid, we again will have a stack of parallel contours approximating the contour 
surface, each planar set of contours corresponding to a particular y coordinate. 
Likewise, if we contour all y-z slices of the three-dimensional grid, we will have a 
stack of parallel contours approximating the contour surface, each planar set of 
contours corresponding to a particular x coordinate. The assemblage of the 
three sets of parallel, planar contours, i.e. the simultaneous display of all the 
contours created for the x-y, x-z, and y-z planes of the three-dimensional grid, 
produces a "chicken-wire-like’' contour surface display (see Figure 1). The 
three-dimensional contour surface display described in this study is created by 
such a procedure. 

Given that the core of the contour surface display generation algorithm is 
the two-dimensional slice of the three-dimensional grid, it is best that we start 
our study with an understanding of the operations performed on that slice. Fig- 
ure 2 shows a single, x-y, two-dimensional grid, with the contours drawn 
corresponding to contour level 50. Figure 3 shows that same two-dimensional 
grid, with the contours drawn corresponding to contour level 100. The two- 
dimensional grid of those figures is a 4 x 5 grid; it has four values in the x direc- 
tion, and five values in the y direction. The goal of the two-dimensional contour- 
ing operation for such a grid is the determination of where lines are drawn on 
that grid given a fixed contour level k. In order to develop an intuitive feel for 
that determination mechanism, we restrict our focus to a small portion of the 
complete two-dimensional grid, the 2x2 subgrid. The 2x2 subgrid is defined to 
be that portion of the two-dimensional grid bounded by four adjacent grid 
points. In the two-dimensional grid of Figures 2 and 3, the lower, lefthand 2x2 
subgrid is bounded by points (1,1), (2,1), (2,2), and (1,2). The upper righthand 
2x2 subgrid of the same example is bounded by points (3,4), (4,4), (4,5), and 
(3.5). 

2. A Model for Contouring the 2x2 Subgrid 

The procedure used to generate the contours for a single 2x2 subgrid is 
the core part of two-dimensional contouring. If we compute the contours 
corresponding to contour level k for all 2 x 2 subgrids of a two-dimensional grid, 
then we will have determined the complete set of contours for that grid Note 




Figure 2 

Example Contour Grid with Contours Drawn for Level 50 



va 



ia 



40 



10 ..... / .. 


193 \ 


50 __ 


3 


30 


30 


50 / 


150. \ 


40 . _ 


20 


5Q 


is 



13 



23 



30 



70 



Figure 3 

Example Contour Grid with Contours Drawn for Level 100 
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that this does not make any statement as to the efficiency of that picture, i.e. 
there can be duplicate copies of contours, particularly for contours drawn along 
the border of a 2 x 2 subgrid (shared edge). In order to provide an intuitive feel 
for contour generation on the 2x2 subgrid, we briefly summarize the operations 
that comprise that procedure in order to highlight potential problems. 

The procedure used to generate the contours for a particular 2x2 subgrid 
first determines if any contours should be generated for that subgrid. That 
determination is based upon whether any of the subgrid’s edges contain the 
desired contour level k. An edge contains contour level k if the value of that 
contour level is within the range of values defined by the grid points that 
comprise the edge. 

The next part of the contour generation procedure for the 2x2 subgrid is 
the computation of the contour edge intersections for any subgrid edges shown 
to contain the contour level. The point of intersection is computed through 
linear interpolation, using the grid values assigned to the endpoints of the edge 
and their corresponding coordinates. The point of intersection represents the 
location on the subgrid edge corresponding to the contour level k. 

The determination of the connectivity necessary to form the appropriate 
contours from the list of edge intersections is the next part of the contour gen- 
eration procedure. Before attempting to describe the procedure that assigns 
those connectivities, we first examine the subgrid’s contour crossing possibili- 
ties. We accomplish that by looking at Figure 4, which shows all possible ways 
for contours to cross or intersect a 2 x 2 subgrid. 

In Figure 4, there are ten cases, each of which belongs to one of three con- 
tour crossing categories: (1) single edge crossings of the 2x2, (2) double edge 
crossings of the 2x2, and (3) constant edge borders at the contour level for the 
2x2. The ten cases are drawn according to the following small set of rules for 
contour crossings. 

(1) Contours are directed by the values associated with the edges, and 
are directed towards edge intersections. 

(2) For non-equivalued edges, if contours are indicated for a particular 
2x2 subgrid, i.e. there are edges in the subgrid that contain the con- 
tour level, there is only one point of intersection for each edge of 
that subgrid. 

(3) Contours are continuous, Le. if a contour enters a 2 x 2 subgrid, it 
must also leave that 2x2 subgrid. 

(4) Equivalued subgrid edges at the contour level are special cases, and 
are drawn in their entirety. The only exception to this rule is that con- 
stant valued 2x2 subgrids are not drawn. This is by convention. 

The first rule, that of contours being directed by the values associated with 
the edges, and contours being directed towards edge intersections, means that 
one determines the placement of contours, and hence, the connectivity of the 
edge intersections, by using both the values assigned to the endpoints of each 
edge of the subgrid, and the computed intersections of that subgrid. The impor- 
tance of this rule is twofold. First, it indicates that no outside forces or parame- 
ters direct contour placement. Second, the rule indicates that computed inter- 
sections are not the sole basis for determining the connectivity of the contours. 
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The second rule, that of there being only one point of intersection for each 
edge of the 2x2 subgrid, means that for an edge intersected by a contour, 
there are no other points along that edge with contour value k. Since the range 
of values assigned to that edge is continuous and monotonic, the origin of this 
rule is clear. Note that this rule does not apply to equivalued edges. 

The third rule, that of contours being continuous, means that if a contour 
enters a 2 x 2 subgrid, it must also leave that subgrid, i.e. the contour does not 
terminate inexplicably in the middle of the 2x2. A corollary of this rule, in 
combination with the first and second rules, states that contours entering a 
subgrid through an edge must leave via a different edge. Again, note that this 
corollary does not apply to cases with equivalued edges. This corollary holds 
though for contours tangent to grid points of the subgrid if we arbitrarily assign 
one of the subgrid edges to that grid point as the entering edge, and assign the 
other edge as the leaving edge. 

The last rule, that of equivalued subgrid edges at the contour level being 
special cases and being drawn in their entirety, is a rule based on visual expec- 
tations for the contour for such a subgrid edge. The only exception to this rule, 
that of not drawing constant valued 2x2 subgrids, is adopted by general con- 
vention. 

Once we have an idea of the types of contour crossings possible for a 2 x 2 
subgrid, and once we have an outline of the rules used in composing those possi- 
bilities, we can then address the problem of forming a procedure for assigning 
connectivities to the computed edge intersections. Starting with the simplest 
cases of Figure 4, the equivalued edge cases, we clearly see that the connectivity 
generation procedure for subgrids containing such edges at the contour level is 
relatively simple once those equivalued edges have been detected. If we find 
that we have a "constant 2 x 2", we do not need to issue any coordinates or con- 
nectivities because by convention we have decided not to draw that case. The 
other two possibilities, the "contour along one edge", or the "contours along two 
edges" cases, are equally as simple. The only operation necessary once such 
cases have been detected is to issue coordinates and connectivities correspond- 
ing to the detected edges. 

At first glance, given the edge intersections for a 2 x 2 subgrid, the connec- 
tivity generation procedure for the single contour cases of Figure 4 seems quite 
easy. It appears as if the only operation that has to be done is to issue coordi- 
nates and connectivities corresponding to the straight line between the two 
points of edge inter sec tion. Such a procedure works well if we know that we 
have a single contour crossing the subgrid. The only single contour crossing 
case for which this will not work is the "contour tangent to the 2 x 2" case, which 
is an even simpler case for connectivity generation. 

It is not until we consider the two contours crossing the subgrid cases of 
Figure 4 that we realize the potential for problems with the above single contour 
crossing procedure. A procedure based only on connecting edge intersections 
cannot differentiate between cases such as the "two contours tangent to the 
2 x 2", and the "contour across the diagonal" cases. There are other similar con- 
nectivity generation problems evident for the two contours crossing cases. The 
"two contours through adjacent edges" case has four edge intersections. For 
that case, information needs to be provided to the connectivity generation pro- 
cedure that determines which of three possible intersection pairs should be con- 
nected. 

Now that we have established a background for the connectivity problem for 
contour crossings of the 2x2 subgrid, we can detail the procedure used to 
solve that problem. Before we describe that algorithm, we first briefly review 
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some of the problems cited Ln the literature for two-dimensional contouring. 

3. Two-Dimensional Contouring Literature Problems 

The literature on two-dimensional contouring, and the use of two- 
dimensional contouring for creating a contour surface display is extensive, 
encompassing a number of fields ([Zyda, 1984], [Zyda, 1983], [Zyda, 1982], 
[Zyda.1981], [Dutton. 1977], [Gold, 1977], [Wright, 1979], [Cottafava, 1969], 
[Dayhoff, 1963], [Faber, 1979], [McLain, 1974], [Sabin. 1980], [Sutcliffe, 1980], 
[Sutcliffe, 1976a], and [Sutcliffe, 1976b]. A thorough review of the historical 
development of two-dimensional contouring algorithms and the properties of 
those algorithms is found in [Sutcliffe, 1980]. Many of the contouring algorithms 
presented in that study are flawed either in that they generate an incorrect pic- 
ture for some contour crossing cases, or in that they require special handling 
for "problem” 2x2 subgrids. Some of the typical algorithm problems detailed 
are identical to those described above, i.e. they concern “degenerate points”, 
grid points that are equal to the contour level, or "saddle points", 2x2 subgrids 
where there are ambiguities as to which points to connect. In all of the algo- 
rithms reviewed, no attempt is made to fit the special cases inside of a general 
algorithmic framework. This is quite evident for the subgrid having a saddle 
point. That contour crossing case is handled by selecting the two lines "for 
which the direction changes the least" when compared against neighboring 
subgrids [Sutcliffe, 1980]. Again, this requires special algorithmic resolution. 
None of the papers attempts to build a general framework useful for the genera- 
tion of the coordinates and drawing instructions for any 2x2 subgrid. The fol- 
lowing section describes both a data structure, the contouring tree, and an algo- 
rithm for using that data structure, that provide both a coherent framework for 
2x2 subgrid contouring, and a comprehensive resolution to the 2x2 subgrid 
crossing problem* 

4. The Contouring Thee 

A contouring tree is a data structure that represents the edge value rela- 
tionships of a 2 x 2 subgrid in a form that permits the rapid generation of the 
contour display for any contour level contained within the represented subgrid 
(see Figure 5). The formulation of the contouring tree is based upon the obser- 
vation that for any two-dimensional grid a continuous series of contour displays 
can be created for contour levels in the range of the minimum and maximum 
grid values (see Figure 6, and [Zyda, 1984], [Zyda 1983], [Zyda 1982], and 
[Zyda 1981]). 

The use of the contouring tree is outlined best with an example of a small 
two-dimensional grid. Figures 2 and 3 depict the contours generated for con- 
tour levels 50 and 100. The contours at level 100 are closed contours, forming 
simple, connected loops. The contours at level 50 are open contours. Figures 5 
and 7 present the contouring trees created for two 2x2 subgrids of the 4x5 
plane. The edges of the contouring trees correspond to the directed, downhill 
edges inscribed on the 2x2 subgrids of the figures. There are eight directed 
edges on each subgrid, four for the boundary edges and four for the edges to the 
subgrid’s center point. The value used for the center point is the average of the 
four values comprising the corners of the 2x2 subgrid. (A reference as to the 
usefulness of the center point average value in generating smooth contours is 
found in [Sutcliffe, 1980].) The edges of the contouring trees are ordered, main- 
taining the same counterclockwise ordering as in the original subgrids. A "1" 
under a node indicates that a setpoint display command should be generated for 
any coordinate that is created along an edge that has that connectivity on its 
lower valued node. A ”0" indicates a drawto display command in a similar 
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Coordinates Generated for Sample 2x2 Subgrid 
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Figure 6 

Example Contour Grid with Contours Drawn for Multiple Contour Levels 
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Coordinates Generated for Sample 2x2 Subgrid with Saddle Point 
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fashion and a "2" indicates a drawpoint. 

Display generation from a contouring tree is accomplished by performing a 
pre-order traversal of that contouring tree, producing a coordinate and drawing 
instruction whenever the desired contour level is found to be within the range of 
an edge of the contouring tree. A pre-order traversal visits the root, the left 
subtree, the middle subtree, and then the right subtree. An edge’s range is 
defined to be the set of values between those associated with the nodes on either 
end of the edge. More precisely, we say a contour level is within an edge if the 
following condition holds: 

lower _jiode’s_yalue <_contour_level < higher_node’s_yalue 

For example, in Figure 5a at contour level 100, we issue coordinates and drawing 
instructions for the edges (2,2)-(3,2), (2,2)-(2.5,2.5), and (2,2)-(2,3). The drawing 
instruction issued for each of these edges is again the one associated with the 
lower valued node of the edge. The coordinate for each of these edges is gen- 
erated by a linear interpolation of the edge’s endpoint coordinates according to 
the decrease in contour level along the edge. The coordinates and drawing 
instructions generated for the contouring trees of Figures 5a and 7a are 
represented in Figures 5b and 7b. 

There are some subtleties not evident from the above that are best detailed 
using a pseudocode description of the traversal algorithm. Figure 8 depicts the 
traversal procedure for the contouring tree assuming a particular data organi- 
zation. The notation is quite standard. The pointers to the descendent nodes of 
NODE are LEFT(NODE), MIDDLE(NODE). and RIGKT(NODE). For each node of the 
contouring tree, there are three pieces of information: the value associated with 
the node, VALUE(NODE), the coordinate associated with the node, XYZ(NODE), 
and the connectivity associated with the node, CONN(NODE). 

The generation of coordinates and drawing instructions from a contouring 
tree begins with routine C0NT0UR_SUBGRID of Figure 8. That routine receives a 
pointer to the root node of the contouring tree. It then starts the traversal by 
calling routine VISIT with that root node. Routine VISIT checks to see if the edge 
defined by the passed in node and that node’s ancestor, NODE and ANCESTOR, 
contains the contour level. If the edge does contain the contour level, the edge 
intersection coordinate is computed using linear interpolation and issued to the 
display along with the connectivity associated with that node, CONN(NODE). If 
we issue a coordinate and connectivity for a node, we need to check the subtree 
under that node for equivalued edges. If an equivalued edge at the contour level 
is found, a coordinate and drawing instruction pair are issued for that 
equivalued edge (routine VISIT _SUBTREE). Once a coordinate and drawing 
instruction pair have been issued for an edge, and once the subtree beneath 
that edge has been investigated for equivalued edges, further traversal of that 
subtree is terminated. If an edge is found not to contain the contour level, the 
traversal continues as depicted at the bottom of routine VISIT. 

The pre-order traversal procedure described above generates the coordi- 
nates and drawing instructions for the part of the 2x2 subgrid the contouring 
tree represents. To generate the coordinates for a larger two-dimensional grid, 
we generate the contouring trees for each 2x2 subgrid of that grid, and then 
apply the traversed procedure to those trees. We note here that no ordering is 
required in the generation of coordinates for the 2x2 subgrids. The coordinate 
and drawing instruction set generated for each 2x2 subgrid is complete and 
independent of the picture generated for any neighboring 2x2. 



Contouring Tree Description 

Pointers to descendent nodes: 

LEFT(NODE) 

MIDDLE(NODE) 

RIGHT(NODE) 



Values associated with each node: 

VALUE(NODE): grid value 

XYZ(NODE) : coordinate of that grid value. 

CONN(NODE) : drawing instruction. 



Procedure CONTOUR_SUBCRID(ROOT) 

V1SIT(R00T,R00T) # begin the traversal of the pointed at 
jf contouring tree. 



end. 



Procedure VIS1T(N0DE, ANCESTOR) 

if{NODE == NULL) 

\ 

return 

j 

if((VALUE(NODE) <= CONTOURJ.EVEL < VALUE(ANCESTOR)) 

OR 

(VALUE(NODE) == CONTOUR_LEVEL AND NODE == ANCESTOR)) 

\ 

§ Edge contains the contour level. 

Issue a coordinate computed via linear interpolation 
along the edge. 

Issue CONN(NODE) as the drawing instruction. 



Figure 8 

Pseudocode of the Traversal Algorithm for the Contouring Tree 



ft Check subtrees of this node for equivalued edges. 
VISIT _SUBTREE(LEFT(NODE),NODE) 
V1SIT_SUBTREE(MIDDLE(N0DE),N0DE) 

VISIT _SUBTREE(RIGHT(NODE),NODE) 

return ft no need to examine the subtree further. 

} ft endif coordinates were generated for an edge. 



VISIT(LEFT(NODE).NODE) ft visit left subtree. 
VISIT(MIDDLE(NODE).NODE) ft visit middle subtree. 
V1SIT(RIGHT(N0DE),N0DE) ft visit right subtree. 

return 

end 



Procedure VISIT_SUBTREE(SUBNODE,SUBANCESTOR) 
if(SUBNODE == NULL) 
return 

\ 

if(VALUE(SUBNODE) == CONTOUR_LEVEL) 

1 



Issue coordinates for the equivalued edge. 
Setpoint on XYZ(SUBANCESTOR). 

Drawto XYZ(SUBNODE). 



I 



V1SIT_SUBTREE(LEFT(SUBN0DE), SUBNODE) 
VISlT_SUBTREE(MIDDLE(SUBNODE). SUBNODE) 
VISIT_SUBTREE(RICHT(SUBNODE), SUBNODE) 

return 



end 



Figure 8 (continued) 

Pseudocode of the Traversal Algorithm for the Contouring Tree 
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4. 1. Contouring Tree Use Discussion 

Having presented the use of the contouring tree, we must look back and 
discuss its capabilities and limitations. The initial impression is that the con- 
touring tree provides a nice, uniform framework for generating the coordinates 
and drawing instructions appropriate to the 2x2 subgrid. The algorithm takes 
care of the single contour crossing cases quite readily. The algorithm also takes 
care of the difficult two contours crossing case for the 2x2 subgrid. The algo- 
rithm correctly handles subgrids containing equivalued lines at the contour 
level. The algorithm also handles subgrids containing a single grid point at the 
contour level. 

The core problems with this algorithm all concern issues of picture 
efficiency. Since the display generated for each 2x2 subgrid is generated 
independently of any neighboring 2x2 subgrids, equivalued lines at the contour 
level on the border of a subgrid will be duplicated. A similar problem occurs for 
subgrid corner values that equal the contour level. If we display either of the 
above cases on a calligraphic display device, we will see a bright line for the 
equivalued edge, and a bright point for the grid value equal to the contour level. 
Another problem, also due to the independent computation of each 2x2 
subgrid, is that no ordering is provided for coordinates that come out of this 
algorithm. For calligraphic displays, this is a problem because for such devices 
electron beam movement is expensive. A contour display that causes the max- 
imum movement of the electron beam every other subgrid greatly decreases 
the the vector capability of the calligraphic display device. 

There are three possible solutions to the first problem, that of duplicate 
vectors. The easiest solution is to choose an output display device for which 
such picture inefficiencies do not matter, i.e. a raster display. Vector ordering 
is also eliminated as a problem with this solution. The second solution to the 
vector duplication problem is to set aside points and lines at the contour level 
that correspond to subgrid boundaries. A final pass at the end of the computa- 
tion for a complete two-dimensional plane could readily cull the duplicates. This 
second solution does nothing for the vector ordering problem If this solution to 
the vector duplicates problem is utilized, one either doesn’t worry about the 
vector ordering, or performs a sort on the vectors. The third solution, and the 
most expensive of the three, is to merge the set of trees generated for the two- 
dimensional grid such that duplicate edges in separate trees are eliminated. 
This solution has the added benefit that the resultant contours are generated in 
an order that solves the beam movement problem. This solution is not 
described in detail here and the reader is referred to [Zyda, 1981] for further 
detail. For this study, the first and simplest solution is assumed, and conse- 
quently, the expected output display device is the raster display. 

5. Contouring Tree Construction 

Contouring tree construction is best understood if we describe that pro- 
cedure in graph theoretic terms. We begin by assuming that we have a graph of 
five nodes, each node being one of the subgrid definition nodes or the center 
node of average value. The eight edges on that graph are the subgrid boundary 
edges, and the edges from each subgrid definition point to the center point of 
average value. We can readily assign directions to each edge of this graph using 
the values assigned to each node. Equivalued edges can be assigned an arbi- 
trary direction. (One such assignment is to make equivalued edges along the 
border point in a counterclockwise fashion, and equivalued edges from the 
center point in towards the center.) With these assumptions, each 2x2 subgrid 
is perceived as a directed graph. The question then becomes, how do we obtain 



it Check subtrees of this node for equivalued edges. 
VISIT_SUBTREE(LEFT(NODE),NODE) 
VISiT_SUBTREE(MIDDLE(NODE).NODE) 
VISIT_SUBTREE(R1GHT(N0DE),N0DE) 

return § no need to examine the subtree further. 

j ff endif coordinates were generated for an edge. 



V1SIT(LEFT(N0DE),N0DE) ff visit left subtree. 
V1SIT(MIDDLE(N0DE), NODE) ff visit middle subtree. 
V1SIT(RIGHT(N0DE).N0DE) # visit right subtree. 

return 



end 



Procedure VISIT_SUBTREE(SUBNODE.SUBANCESTOR) 

if(SUBNODE == NULL) 

\ 

return 

i 

if(VALUE(SU BNODE) == CONTOUR _LEVEL) 

\ 

Issue coordinates for the equivalued edge. 
Setpoint on XYZ(SUBANCESTOR). 

Drawto XYZ(SUBNODE). 



VISIT_SUBTREE(LEFr(SUBNODE). SUBNODE) 
VISIT_SUBTREE(MIDDLE(SUBNODE).SUENODE) 
VISIT_SUBTREE(RIGHT(SUBNODE). SUBNODE) 

return 



end 



Figure B (continued) 

Pseudocode of the Traversal Algorithm for the Contouring Tree 
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4. 1. Contouring Tree Use Discussion 

Having presented the use of the contouring tree, we must Look back and 
discuss its capabilities and limitations. The initial impression is that the con- 
touring tree provides a nice, uniform framework for generating the coordinates 
and drawing instructions appropriate to the 2x2 subgrid. The algorithm takes 
care of the single contour crossing cases quite readily. The algorithm also takes 
care of the difficult two contours crossing case for the 2x2 subgrid. The algo- 
rithm correctly handles subgrids containing equivalued lines at the contour 
level. The algorithm also handles subgrids containing a single grid point at the 
contour level. 

The core problems with this algorithm all concern issues of picture 
efficiency. Since the display generated for each 2x2 subgrid is generated 
independently of any neighboring 2x2 subgrids, equivalued lines at the contour 
level on the border of a subgrid will be duplicated. A similar problem occurs for 
subgrid corner values that equal the contour level. If we display either of the 
above cases on a calligraphic display device, we will see a bright line for the 
equivalued edge, and a bright point for the grid value equal to the contour level. 
Another problem, also due to the independent computation of each 2x2 
subgrid, is that no ordering is provided for coordinates that come out of this 
algorithm. For calligraphic displays, this is a problem because for such devices 
electron beam movement is expensive. A contour display that causes the max- 
imum movement of the electron beam every other subgrid greatly decreases 
the the vector capability of the calligraphic display device. 

There are three possible solutions to the first problem, that of duplicate 
vectors. The easiest solution is to choose an output display device for which 
such picture inefficiencies do not matter, i.e. a raster display. Vector ordering 
is also eliminated as a problem with this solution. The second solution to the 
vector duplication problem is to set aside points and lines at the contour level 
that correspond to subgrid boundaries. A final pass at the end of the computa- 
tion for a complete two-dimensional plane could readily cull the duplicates. This 
second solution does nothing for the vector ordering problem If this solution to 
the vector duplicates problem is utilized, one either doesn’t worry about the 
vector ordering, or performs a sort on the vectors. The third solution, and the 
most expensive of the three, is to merge the set of trees generated for the two- 
dimensional grid such that duplicate edges in separate trees are eliminated. 
This solution has the added benefit that the resultant contours are generated in 
an order that solves the beam movement problem. This solution is not 
described in detail here and the reader is referred to [Zyda, 1981] for further 
detail. For this study, the first and simplest solution is assumed, and conse- 
quently, the expected output display device is the raster display. 

5. Contouring Tree Construction 

Contouring tree construction is best understood if we describe that pro- 
cedure in graph theoretic terms. We begin by assuming that we have a graph of 
five nodes, each node being one of the subgrid definition nodes or the center 
node of average value. The eight edges on that graph are the subgrid boundary 
edges, and the edges from each subgrid definition point to the center point of 
average value. We can readily assign directions to each edge of this graph using 
the values assigned to each node. Equivalued edges can be assigned an arbi- 
trary direction. (One such assignment is to make equivalued edges along the 
border point in a counterclockwise fashion, and equivalued edges from the 
center point in towards the center.) With these assumptions, each 2x2 subgrid 
is perceived as a directed graph. The question then becomes, how do we obtain 
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the contouring tree, or trees from this directed graph? We can put this question 
In terms of graph theory if we notice that a contouring tree is a directed tree. 
The problem then becomes one of obtaining the directed tree, or trees, from the 
directed graph such that the order of edge attachment in the tree corresponds 
to the order in the directed graph (the 2x2 subgrid). From graph theory, we 
have the requirement that a directed tree has the in-degree of its root node 
equal to zero, and the in-degree of every other node equal to one [Even, 1979]. 
To examine the in-degree of each node of the directed graph, we must construct 
the in-degree matrix D for that graph. The in-degree matrix D of a directed 
graph G is defined in [Even, 1979] as: 

D(i,j) = in-degree (i), if i = j. 

D(i.j) — -k, if i is not equal to j, 

where k is the number of 
edges in G from i to j 
(i.e., -1 for all our graphs). 

Figures 9 and 10 show the in-degree matricies for our example 2x2 subgrids. 

From the figures, we note that the roots of the contouring trees are recog- 
nizable from D as D(i.i) = 0. This matches the first part of the directed tree 
requirement. Further examination of the diagonal of the in-degree matricies 
introduces two difficulties in our attempts to convert the directed graphs into 
directed trees: (l) multiple roots (in-degree(v) = 0 for more than one node) and 
(2) in-degree(v) > 1 for some nodes v. (Note: we have assumed that we have the 
structure represented by the directed graph and that we can manipulate it.) 

The first problem, that of multiple roots, is handled by producing multiple 
sets of verticies and multiple in-degree matricies such that there is only one 
root per in-degree matrix. For our case, the maximum number of roots for a 
single 2x2 subgrid is two. We eliminate this problem by alternately removing 
each root node from the complete set of verticies and all edges attached to that 
deleted node and then making separate in-degree matricies. The second prob- 
lem in the conversion of the directed graph into a directed tree, that of in- 
degree^) > 1 for some nodes, is resolved by node duplication. For each diago- 
nal entry D(i,i) = n, where n > 1, we create n-1 duplicates of that node, for a 
total of n, taking care to copy the appropriate values, coordinates, etc. We then 
reassign the original edges that went to the single node, such that each edge 
receives its own copy of the duplicated node. The edges that are reassigned are 
those between each node of column i of D that has a -1 and each of the n dupli- 
cate nodes. When performed for each in-degree matrix created, this operation 
forms a new directed graph that is the desired directed tree 

5.1. Drawing Command Placement 

The above section detailed the creation of the contouring trees for the 2x2 
subgrid utilizing nothing more than basic graph definitions, and simple con- 
struction procedures. The only thing remaining for completing this construc- 
tion is the placement of the drawing commands into the contouring trees. 

Drawing commands are placed in the contouring tree to indicate when a 
line enters the region represented by the contouring tree either from a neigh- 
boring subgrid or from a location off of the grid. If we look at the structure of 
the contouring tree, such as that exhibited in Figure 5a, and consider that dur- 
ing the traversal, the edges are examined in a counterclockwise, and downward 
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Figure 9 

In-degree Matrix for the Directed Graph Superimposed on a 2 x 2 Subgrid 
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Figure 10 

In-degree Matrix for the Directed Graph Superimposed 
on a 2 x 2 Subgrid with Saddle Point 
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ordering from the root, we note that we need to place setpoint drawing com- 
mands on the lower valued node of each edge that presents a new lowest value 
for the tree. (Note that the drawing command setpoint indicates to the display 
device that it should move its "drawing instrument”, i.e. electron beam, pen, 
etc., in a non-drawing mode to the specified location, and that it should then 
place that drawing instrument into a drawing mode. Drawto Indicates to the 
display device that it should move its drawing instrument in a drawing mode to 
the specified location. Drawpoint indicates to the display device that it should 
move its drawing instrument in a non-drawing mode to the specified location, 
and that it should then turn that drawing instrument on for the space of a single 
point.) We insert these drawing commands by way of a pre-order traversal of 
the directed tree, placing a setpoint command on each node that is a new lowest 
value for the tree. This drawing command placement strategy is based upon the 
fact that if we have a contour level for which we desire a picture, the first draw- 
ing command we generate for any contouring tree is a setpoint. Although fairly 
effective, this procedure does not provide a complete solution to drawing com- 
mand insertion. Some neighboring edges in the contouring tree, i.e. edges shar- 
ing an ancestor node, have a "split" between them, i.e., the edges axe not 
immediate counterclockwise neighbors in the original grid. In this case, we 
must indicate the discontinuity in the contouring tree. We register the discon- 
tinuity on the lower valued node of the edge where the discontinuity occurs. For 
example, in Figure 5a the edges (3,3)-{3,2) and (3,3)-(2,3) are neighbors in the 
contouring tree but are not immediate neighbors in the original grid. We indi- 
cate this split by placing a "1" on the lower valued node of edge (3,3)-{2,3). 

In order to recognize the nodes that require a drawing command indicating 
a split edge in the contouring tree, we must first examine where split edges 
occur in the 2x2 subgrid. These occur in the 2x2 subgrid wherever the 
subgrid has edge directions and tree edge configurations as depicted in Figure 
11, Le. where there are neighboring tree edges not directly corresponding to 
neighboring subgrid edges. There are two cases in the figure. 

In case 1, the split edge is edge (B,D). This edge neighbors edge (B.C) in the 
contouring tree but is not the immediate counterclockwise neighbor edge of 
(B,C). The lower valued node of edge (B.D), D, receives a setpoint drawing com- 
mand. From the figure, we see that node D has a possibility of being a "sink", 
Le., a node with all incoming edges. This depends upon the the direction of edge 
(A,D). There are two cases to consider (la) A — > D or (lb) D — > A. Case la, 
A — > D, is easy to show as possible because it may be seen in Figure 5a. Case 
lb, D — > A, is somewhat more difficult, because we must show that it cannot 
possibly occur within the relations specified on Figure 11. We begin by compiling 
a small set of the given relations on the directed graph shown in the figure: 



R £_C 
R > A 
R> D 
B £_D 
A £_D 



Case lb then becomes the question, is it possible to get D >.A? Assume that 
D >.A is possible. Now we have B £_D so we can replace D with B. We get B >_A, 
which is a contradiction of our original given, A £_B, unless A = B. So we cannot 
have D — > A unless A = B. Can we eliminate the possibility of the edge from 
A to B pointing from the center outward in the case where A = B? We can do this 
when we set up the original directed graph, biasing edge selection of constant 



CASE 1 



0 B 




R C 




CASE 2 






C D 



ALL SPLIT EDGE POSSIBILITIES FOR THE 2 X 2 SUBGRID 
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valued edges to point always towards the center. Since the direction assigned to 
a constant edge is, in fact, arbitrary, we can assume that we never see the con- 
dition D — > A, given the set of fixed directions assigned to case 1. Hence, for 
case 1 of Figure 11, a split edge node is indicated during pre-order traversal of 
the directed tree whenever we encounter, for the first time, a node represent- 
ing a sink grid point. We now show a similar finding for case 2. 

In case 2 of Figure 11, the split edge is edge (B.D). This edge neighbors 
edge (B.C) in the contouring tree but is not the immediate counterclockwise 
neighbor edge of (B.C). The lower valued node of edge (B.D), D, receives a set- 
point drawing command. From the figure, we see that node D is possibly a sink. 
This depends on the direction of edge (A,D). There are two cases to consider: 
(2a) A — > D or (2b) D — > A. Case 2a, A — > D, occurs in Figure 7a for the two 
maxima example so we know this case is possible. Case 2b, D — > A, is more 
difficult because we must show that it cannot occur. We list a few relations 
describing the partial directed graph of the figure: 

R£_C 
R>B 
R> D 
B > C 
A> B 
B > D 



The question that arises is whether it is possible to get D >.A? Assume D >_A is 
possible. Now we have B >_D so we can replace D with B. This gives B >.A, which 
is a contradiction of our original given, A >_B, unless A = B. So we cannot have 
D — > A unless A = B. Clearly we already have eliminated the possibility of the 
edge pointing from B to A by our initial configuration and by our biased edge 
selection for constant edges. Hence, for case 2 of Figure 11, a split edge node is 
indicated during pre-order traversal of the directed tree whenever we 
encounter, for the first time, a node representing a sink grid point. 

There are no other split edge configurations possible in the 2x2 subgrid. If 
we add a procedure that checks not only the new lowest value in the directed 
tree but also the first encounter with a sink grid node during the same pre-order 
traversal, we will correctly place the drawing commands in the directed tree, 
converting it into the desired contouring tree. 

6. Completeness for the Subgrid Contouring Algorithm 

An algorithm is complete with respect to the subgrid contour generation 
problem if it always generates the expected picture for any proffered subgrid. 
Since it is impossible to test a subgrid contouring algorithm by trying it against 
all possible subgrids, we must rely on another mechanism. Our ability to 
characterize the total number of possible contour crossing cases in as small a 
number as ten brings to mind such a mechanism for validating the contouring 
tree algorithm. 



For each 2x2 case possibility, show that all possible contouring trees 
generatable from the relations specified for the case either 



(1) produce the equivalent picture (to the 
expected), 

(2) or are in actuality the contouring trees of 
another case, i.e. we show that the equivalent set of 
coordinates is impossible to generate from the con- 
touring tree. 

For each subgrid crossing case... 

i 



For each boundary point of the 2x2 subgrid... 

\ 



Use the boundary point as the root of the tree. 

Generate the 32 possible trees that have that root (and ini- 
tial tree configuration). 

For each tree... 

\ 



Can we obtain a set of coordinates and drawing 
instructions equivalent to the expected set for this 
2x2 case? 

If we have the equivalent picture, then we are OK. Next. 

If we can’t get the equivalent picture, we must show 
that we violate the relations set up for the problem. 
This implies a different case. Next. If we don’t violate 
the relations, that implies a problem with the contour- 
ing tree algorithm. 




i 



Figure 12 

Outline of the Contouring Tree Algorithm Completeness Proof 
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If we are able to show for each of the ten cases of Figure 4 that all possible 
contouring trees generatable from the relations specified for each case either 

(1) produce a picture equivalent to that expected for the 
case, 

(2) or are in actuality the contouring trees of a different case, 
i.e. we show that the expected set of coordinates and drawing 
instructions is impossible to generate from the contouring 
tree (we violate the given relations), 

then we will have shown completeness for the algorithm. This is a fairly large 
proof, an outline of which can be seen in Figure 12. 

The idea behind the proof for an individual case is to use each boundary 
point of the subgrid as the root of the set of all possible trees generatable from 
that root. When we fix the root for a contouring tree, we fix the three immediate 
descendents of that root and the three edges to those descendents. This leaves 
unspecified the placement of five edges in the contouring tree. If we assume 
that each edge direction is possible, this means that once we fix the root of a 
tree, there are 32 possible trees that can be generated. If we cycle through 
each boundary point as the root of the tree, there are 32 times 4 possible trees 
for each case. For each of the 129 possible trees, we then check to see if we can 
obtain a set of coordinates and drawing instructions equivalent to the expected 
set for the case (see Figure 4). If we obtain the equivalent picture from a tree, 
we go on to the next tree (or case). If we cannot obtain the equivalent picture, 
we must show that we violate the relations given for the case. Such a violation 
implies that we have a different case, and hence, can go on to the next tree (or 
case). The negative indications for the proof are if we cannot get the equivalent 
picture but do not violate the relations set up for the case. Such an indication 
would imply a problem with the contouring tree algorithm. 

The magnitude of this proof is immense. There are ten different subgrid 
cases, each with 128 possible contouring trees. The evaluation of this total set of 
trees has been performed utilizing a program that systematically traverses each 
possible tree, checking for either the equivalent picture, or a violation of the 
given relations. The results of this study indicate that the contouring algorithm 
is complete. 

7. Conclusions 

This study has described an algorithm for generating the contour lines and 
drawing instructions for a single 2x2 subgrid in a manner that is independent 
from the calculations required for any other 2x2 subgrid. This algorithm has 
been shown to be part of a larger algorithm used to generate the contours for 
two-dimensioneil grids composed of multiple 2x2 subgrids. It has edso been 
shown to be part of an even larger algorithm used to generate the contour sur- 
face display for three-dimensional grids composed of multiple, parallel two- 
dimensional grids. The fact that we can generate the contour lines for each 
2x2 subgrid independently from the calculations required for any other 2x2 
subgrid means that our algorithm for the generation of contour surface displays 
is highly decomposable. This becomes quite striking if we consider that a typi- 
cal 30 x 30 x 30 three-dimensional grid has 75,690 2x2 subgrids. For such a 
case, there is a potential for 75,690 concurrent computations. Given a technol- 
ogy that allows such parallelism and assuming that we can both load and unload 
the data from that system, it appears as if this algorithm has the potential for 
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speeding up one of the most frequently used graphics computations. 
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